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METHOD AND APPARATUS FOR PRODUCING AN IMAGE CONTAINING 

DEPTH INFORMATION 

Field of the Invention 
5 This invention relates to a method and apparatus for 
producing an image containing depth information. The 
invention has particular application to thermal images in 
the 8 to 12 micron band of the electromagnetic spectrum. 
However, the method and apparatus is also useful with 
10 images taken in other parts of the electromagnetic 

spectrum (and in particular visible light), as well as 
electron beam images, and other forms of particle 
radiation images as well as acoustic images. 

15 Background Art 

Thermal images are used to produce a visual image of a 
scene in darkness. In such images, articles in the scene 
which are above ambient temperature can be easily 
visualised. However, one major difficulty with 

20 conventional thermal images is that the images contain 
little or no depth information, so that it is very 
difficult to ascertain whether one part of the image is in 
front of or behind another part of the image. Thus, when 
a thermal image is taken of a scene which, for example, 

25 may contain people, the thermal image will clearly show 
each of the individual people, but will not provide 
sufficient information to ascertain the relative positions 
of the people in space. 

3 0 In many applications, it is necessary that the image does 
include depth information so that the image provides 
sufficient information of the scene to suit a particular 
purpose . 

3 5 Summary of the Invention 

The object of the invention is to provide a method and 
apparatus of producing images containing improved depth 
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information so that a determination can be made as to the 
distance of various parts of the image relative to one 
another and from the location where the image is taken. 

The invention may be said to reside in a method of 
producing an image containing depth information, 
comprising: 

detecting a radiation wavef ield emanating from 
the scene and taking at least two images of the scene at 
different planes relative to the scene to produce a first 
image comprised of a first set of intensity data values 
and a second image comprised of a second set of intensity 
data values; 

determining an intensity variation of data values 
in the first set of data values relative to other data 
values in the first set of data values to produce a first 
set of intensity variances, and determining an intensity 
variation of data values in the second set of values 
relative to other data values in the second set of values 
to produce a second set of intensity variance data; 

processing the first and second sets of intensity 
variance data to obtain image data of the scene containing 
depth information; and 

coding the image data which have the same depth 
information with a code reference to identify the 
different depth information in the image data. 

Thus, by producing the variance data and then processing 
the variance data, the variance data is able to result in 
the production of depth information in the image which can 
then be coded so that the different depth information is 
apparent in the image. 

In the preferred embodiment of the invention, the 
intensity variance data is processed in accordance with 
the phase image production method disclosed in our co- 
pending International Patent Application Nos. 
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PCT/AU99/00949 and PCT/AU02/01398 . The contents of these 
two International applications are incorporated into this 
specification by this reference. 

5 OJhus, the step of processing the intensity variance data 
preferably comprises: 

(a) producing a representative measure of the 
rate of change of variance data values of said radiatipn 
wave field over a selected surface extending generally 

10 across the wave field; 

(b) producing a representative measure of said 
radiation wave field relating to the scene over said 
selected surface; 

(c) transforming said measure of rate of change 
15 of variance to produce a first integral transform 

representation and applying to said first integral 
transform representation a first filter corresponding to 
the inversion of a first differential operator reflected 
in said measure of rate of change of variance to produce a 
20 first modified integral transform representation; 

(d) applying an inverse of said first integral 
transform to said first modified integral transform 
representation to produce an untransformed representation; 

(e) applying a correction based on said measure 
25 over said selected surface to said untransformed 

representation ; 

(f ) transforming the corrected untransformed 
representation to produce a second integral transform 
representation and applying to said second integral 

3 0 transform representation a second filter corresponding to 
the inversion of a second differential operator reflected 
in the corrected untransformed representation to produce a 
second modified integral transform representation; 

(g) applying an inverse of said second integral 
3 5 transform to said second modified integral transform 

representation to produce . a measure of phase of said 
radiation wave field across said selected plane so as to 
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produce said image data as phase image data containing the 
depth information. 

The step of producing a representative measure of said 
5 radiation wavef ield relating to the scene over the 

selected surface may use intensity- values to obtain the 
representative measure or variance values. 

In particular, in one embodiment, intensity values at the 
10 selected surface are used. 

In another embodiment, values are taken representing 
maximum focus from any of the intensity data values. 

15 In another embodiment, variance data values are used as a 
representative measure. 

In a still further embodiment, maximum variance values 
taken from the sets of variance data values are used. 

20 ' 

Preferably at least one of the first or second 
differential operator has a form based on an optical 
system used to acquire the radiation for producing the 
representative measure of the rate of change of variance 
25 of the radiation wavef ield over the selected surface 
extending generally across the wavef ield. 

Preferably both the first and second differential 
operators have a form based on the optical system. 

30 • 

Preferably the first and second integral transforms are 
produced using a Fourier transform. 

In one embodiment of the invention the differential 
3 5 operators have the forms 
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T P +oe z 
where, 

T p {p)=2M&jrjT^(p,Tj)dTj 
and , 



10 



20 



25 



30 



2np 



£1* 



and wherein /? is radial lateral spatial frecjuency, the r| 
longitudinal spatial frequency and Sz is the defocus 
distance in the plane of the object. Also 

NA 



AM 

* y ^objective 

where NA condensor and NA objective are respectively the numerical 
aperture of the condenser and the objective (These are 
settings and dimensions on the microscope) . p obj is the 
maximum spatial frequency accepted by the objective* 

Preferably the step of taking at least two images of a 
scene comprises taking the first image at a first 
defocused plane to produce the first image as negatively 
focused image data and taking the second image of the 
scene at a second defocused plane to produce positive 
defocused image data, the negative and positive defocused 
image being taken on respective sides of a focal plane 
which would produce an in focus image of the scene. 

Preferably the step of determining an intensity variance 
comprises determining a measure of the sharpness of each 
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of the data values relative to the sharpness of data 
values surrounding that data value* 



In the preferred embodiment of the invention the images 



intensity value is compared with the intensity of 
surrounding pixels in order to provide a variance value at 
each pixel. 

Preferably the variance is obtained by the following 
equation for each pixel: 



wherein n is the particular pixel, and the values 1 to 9 
represent the pixels which surround that pixel in an array 
of pixels. 

Thus, in this embodiment, each pixel is compared with 
pixels which immediately surround the particular pixel in 
order to determine the variance. In other embodiments, 
each pixel could be compared with the pixels which 
surround that pixel, together with the pixels which 
surround the pixel and are one further pixel away and so 
on. The variance could also be determined by other 
methods so that a variance of each pixel relative to the 
other pixels in each data set is obtained. 

Preferably the step of coding parts of the image data 
comprises applying different colours to parts of the image 
which have the same depth information so as to produce a 
visual image in which the relative distance of parts of 
the scene compared to one another can be determined. 



are captured by a charge coupled device and an intensity 
value of pixels in the image is determined, and that 



var, 
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Preferably a third image of the scene is taken at an in 
focus plane to produce a third set of intensity data 
values and the third set of intensity data values are 
5 overlaid with the coded image data to produce a visual 
image containing structural information of the scene as 
well as the different depth information in the scene. 

However, in other embodiments, rather than providing a 
10 visual image, the image could be an electronic image 

stored in a computer, and the code reference could simply 
be a code value allocated by the computer to the different 
depth information so that if images at a particular 
location are required, information as to those parts of 
15 the image can be extracted by the computer by the code 
reference • 

Preferably the different depth information is provided by 
allocated a grey scale value to pixel values in the image. 

20 

The invention may also be said to reside in a computer 
program for performing the above-mentioned method. 

The invention may also be said to reside in an apparatus 
25 for producing an image containing depth information 
comprising: 

a camera for detecting a radiation wavef ield 
emanating from a scene and taking at least two images of 
the scene at different planes relative to the scene to 

3 0 produce a first set of intensity data values and a second 
set of intensity data values; 

a processor for determining an intensity 
variation of data values in the first set of data values 
compared to other data values in the first set of data 

35 values to produce a first set of variance data values, and 
for determining a variance of one data value in the second 
set of data values compared to other data values in the 
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second set of data values to produce a second set of 
variance data values; 

the processor also being for processing the first 
and second variance data values to produce image data of 
5 the scene containing depth information; and 

the processor also being for coding parts of the 
image which have the same depth information with a code 
reference to identify the different depth information in 
the image. 

Thus, the processor for determining an intensity variation 
is a processor for: 

(a) producing a representative measure of the 
rate of change of variance data values of said radiation 
wave field over a selected surface extending generally 
across the wave field; 

(b) producing a representative measure of said 
radiation wave field relating to the scene over said 
selected surface; 

(c) transforming said measure of rate of change 
of variance to produce a first integral transform 
representation and applying to said first integral 
transform representation a first filter corresponding to 
the inversion of a first differential operator reflected 
in said measure of rate of change of variance to produce a 
first modified integral transform representation; 

(d) applying an inverse of said first integral 
transform to said first modified integral transform 
representation to produce an untransf ormed representation; 

(e) applying a correction based on said measure 
over said selected surface to said untransf ormed 
repr es ent at ion ; 

(f) transforming the corrected untransf ormed 
representation to produce a second integral transform 
representation and applying to said second integral 
transform representation a second filter corresponding to 
the inversion of a second differential operator reflected 
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in the corrected untrans formed representation to produce a 
second modified integral transform representation; 

(g) applying an inverse of said second integral 
transform to said second modified integral transform 
representation to produce a measure of phase of said 
radiation wave field across said selected plane so as to 
produce said image data as phase image data containing the 
depth information. 

The step of producing a representative measure of said 
radiation wavef ield relating to the scene over the 
selected surface may use intensity values to obtain the 
representative measure or variance values. 

In particular , in one embodiment , intensity values at the 
selected surface are used. 

In another embodiment/ values are taken representing 
maximum focus from any of the intensity data values. 

In another embodiment, variance data values are used as a 
representative measure - 

In a still further embodiment, maximum variance values 
taken from the sets of variance data values are used. 

Preferably at least one of the first or second 
differential operator has a form based on an optical 
system used to acquire the radiation for producing the 
representative measure of the rate of change of variance 
of the radiation wavef ield over the selected surface 
extending generally across the wavef ield. 

Preferably both the first and second differential 
operators have a form based on the optical system. 

Preferably the first and second integral transforms are 
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10 



15 



5 



produced using a Fourier transform. 

In one embodiment of the invention the differential 
operators have the form: 



where, 

T P {p)=2mdz lifT™{p,tj)dTi 
and 



Tip,*,)- 



27ZP 



and wherein /> is radial lateral spatial frequency, the n 
longitudinal spatial frequency and & is the defocus 
distance in the plane of the object. Also 



NA. 



NA 



objective 



where NA condensor and NA obJeciive are respectively the numerical 
aperture of the condenser and the objective (These are 
settings and dimensions on a microscope, if a microscope 
is used in the capture of the various images by camera) . 



J obj 



is the maximum spatial frequency accepted by the 



objective. 

If images are captured through a different optical system 
such as a telescope or simply by the optics of a camera, 
then the above operators are modified to suit the various 
optics of that system. 
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Preferably the step of taking at least two images of a 
scene comprises taking a first image at a first defocused 
plane to produce a negatively focused image and taking a 
second image of the scene at a second defocused plane to 
produce a positive defocused image, the negative and 
positive defocused image being taken on respective sides 
of a focal plane which would produce an in focus image of 
the scene.. 

Preferably the step of determining an intensity variance 
comprises determining a measure of the sharpness of each 
of the data values relative to the sharpness of data 
values surrounding that data value • 

In the preferred embodiment of the invention the images 
are captured by a charge coupled device and an intensity 
value of pixels in the image is determined, and that 
intensity value is compared with the intensity of 
surrounding pixels in order to provide a variance at each 
pixel . 

Preferably the variance is obtained by the following 
equation for each pixel: 



1 

var.„ = — x 

. n g 



(/„-/ 2 ) 2 +(/ n -/ 4 ) 2 +(/ n -/ 8 ) 2 



+ i (/ "- /7)2+ W (/ "- /9)2 



wherein n is the particular pixel, and the values 1 to 9 
represent the pixels which surround that pixel in an array 
of pixels and I is the intensity value of the respective 
pixel . 

Preferably the step of coding parts of the processed image 
comprises applying different colours to parts of the image 
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which have the same depth information so as to produce a 
visual image in which the relative distance of parts of 
the scene compared to one another can be determined. 

5 However, in other embodiments, rather than providing a 
visual image, the image could be an electronic image 
stored in a computer, and the code reference could simply 
be a code value allocated by the computer to the different 
depth information so that if images at a particular 
10 location are required, information as to those parts of 

the image can be extracted by the computer by reference to 
the code reference. 



15 



20 



25 



30 



35 



Preferably the different depth information is provided by 
allocated a grey scale value to pixel values in the image. 

Brief D escription of the Drawings 

A preferred embodiment of the invention will be described 
with reference to the accompanying drawings in which: 

Figure 1 is a diagram to assist in illustrating 
the preferred embodiment of the invention; 

Figure 2 is a view of an apparatus embodying the 

invention; 

Figure 3 is a diagram showing pixels of a camera 
used to take the images in the preferred embodiment of the 
invention; 

Figure 4 is a diagram illustrating the production 
of different depth information; 

Figure 5 is a diagram illustrating how different 
depth information is coded according to the preferred 
embodiment of the invention; 

Figure 6 is a flowchart showing an implementation 
of the method of the preferred embodiment; 

Figure 7 is a schematic illustration of an 
arrangement which uses a microscope in the capture of 
images ; 

Figure 8 is a schematic drawing of an exemplary 
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system according to one embodiment of the invention; and 

Figure 9 is a flow chart relating to one 
embodiment of the invention. 

Detailed Description of the Preferred Embodiment 
The preferred embodiment of the invention effectively 
processes images in accordance with the phase 
determination technique disclosed in our International 
Applications PCT/AU99/ 00949 and PCT/AU02/01398 . However, 
before the algorithm disclosed in these applications is 
used to process the data relating to the captured scenes 
to produce the image , the data is first manipulated so as 
to produce a variance data so the variance data is 
processed in accordance with the algorithms rather than 
pure intensity data, as disclosed in the above 
International applications. However, apart from this, the 
method of processing the data is the same as that 
disclosed in the above International applications. 

Figure 1 is a diagram showing an image of a scene at three 
different focal planes. The plane I 0 is an in focus plane 
which will be used to capture an in focus image of the 
scene. The scene may include a number of articles or 
objects such as vehicles, personnel, or the like, and in 
general, the preferred embodiment of the invention is to 
be used for obtaining thermal images. However, the 
invention could be used with other types of images and 
improve the depth information in those images. The 
invention has particular application to thermal images 
because thermal images generally contain little or no 
depth information which enables different parts of the 
images to be compared in terms of distance with other 
parts of the image. Thus, in thermal images, it is very 
difficult to ascertain whether one part of the image is 
closer to or further away from the location from where the 
images were captured, or where parts of the image are 
located in space relative to one another. 
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The plane V. ie a negatively defocused plane within the 
depth of field of the camera 20 (see Figure 2) used to 
capture the images, the plane I 0 is the focal plane, and 
5 the plane V + is a positively defocused plane also within 
the depth of field of the camera 20, but positively 
defocused relative to the plane I 0 . In other words, the 
plane V_ is closer to the scene than the in focus plane I 0 
and the plane V + is further away from the scene than the in 
10 focus plane l 0 . 

Typically, parts of the image taken at the plane V_ which 
happen to be in focus t the plane V. will produce a very 
sharp image point 21 at the plane V_, but will be a very- 
fuzzy defocused image at the plane V + as shown by reference 
22. The waveform which produces the images at 21 and 22 
is shown at 23 in Figure 1. Similarly, part of the image 
24 which is defocused and blurry at the plane V_ may be 
focused at the plane V + as represented by reference 25. 
Once again, the form of the wavefield producing these 
parts of the images is shown at 26. 



15 



20 



As is shown in Figure 2, the camera 20 is used to capture 
images of the scene schematically shown at S in Figure 2 
!5 at the planes V_, I 0 and V + so that three sets of intensity 
data at those three planes is captured by the pixels of 
the camera 20. 

In the preferred embodiment of the invention, the camera 
0 20 is, a cryogenically-cooled, 640 x 480 pixel array 

thermal camera operating in the 8. to 12 micron band with 
adjustable -focus Ge lens. The intensity data captured by 
the camera 20 is processed in a processor 40 and images 
may be displayed oh a display 42. 

5 

The processor 40 first determines a variance value for 
each of the pixels for each of the two images at the 



H.-\Karen\KeeD\Iatia Prov - imaga containing depth information!. doc 27/02/04 



- 16 - 

planes v_ and V + . Thus, two sets of variance data is 
therefore produced. 

Figure 3 shows an example of nine pixels of the pixel 
array forming camera 20 and the variance at the pixel I 5 is 
determined in accordance with the following equation: 



1 

var.< = -x 
3 8 



(/5-/ 2 ) 2 +(/s-/4) 2 + (/ 5 -/ 8 ) 2 
+ (/5 " /6)2+ W (/s " /,)2+ ^ (/5 ~ /3)2 



Thus, a variance value for the fifth pixel in the array of 
Figure 3 is determined. in Figure 3, pixels which are 
arranged at the edge of the array can be set to a value 0 
so that they are not included in the calculations. This 
reduces the size of the image slightly by 2 pixel widths. 
However, it makes calculation much easier because all of 
the variances can be determined in the above mentioned 
manner because all pixels will then be surrounded by 8 
other pixels. Similarly, variance values for each of the 
other pixels in the array is determined for each of the 
image pixels of planes V_ and V + so that two sets of 
variance data is produced. The two sets of variance data 
is then used in the algorithm to produce the phase image 
data of the scene instead of the pure intensity data 
originally captured by the camera 20. The use of the 
variance data provides a measure of the sharpness of the 
image at each of the planes V. and V + and therefore, 
inherently includes a measure of the depth of parts of the 
image relative to other parts of the image. The variance 
data is then processed in the same manner as the intensity 
data, as disclosed in the aforementioned International 
applications, so as to produce a phase image of the scene 
in which depth of field information is included. Because 
the variance data is used to produce the phase image, the 
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image will contain very little detail of the actual 
structure of articles or things in the scene, but will 
contain general outlines showing features of the scene. 
The depth of position information in the processed image 
5 can be given by a grey scale value at each of the pixels 
in the processed image data, which in turn therefore 
provide a relative measure of the part of the scene which 
relates to each pixel relative to other parts of the 
scene, as well as relative to the position of parts of the 
10 scene from where the images were taken. For example, as 
is best shown in Figure 4, something which is very close 
to where the image is taken appears darker, as shown by 
reference 61, and something which is further away appears 
much lighter as at 62* 
15 ■ . 

The grey scale values allocated to each of the pixels in 
the processed phase image can be used to provide an actual 
measure of the distance of the parts of the scene from 
where the image was taken and relative to one another. 
2 0 For example, if the optics of the camera is known and the 
distance between the planes V x and V + is also known, actual 
measurement values can be associated with various grey 
scale values in the image to provide a distance measure of 
parts of the scene relative to one another, and also from 
25 where the images were taken. 

The image shown in Figure 4 can then be coded so as to 
more clearly bring out the different distances of articles 
in the image by applying a code reference to parts of the 

3 0 image which have the same grey scale value. This can be 
done in a variety of different ways, but most preferably 
is performed by coding parts of the image with a false 
colour palette so that parts of the image which have the 
same grey scale value are coloured with the same colour. 

35 This therefore makes it very much easier to perceive depth 
information in the image than in the mere conventional 
thermal image in which the entire image appears to be 
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approximately the same shade. For instance, in Figure 5, 
the grey scale values are shown for near articles in the 
image which may be say 10 metres from where the image is 
taken, up to 100 metres away. Those grey scale values can 
5 be allocated with different colours in the ma nn er referred 
to above , so that a particular colour can be associated 
with a particular distance in the image to make it 
therefore much more easy to perceive where articles in the 
image are actually located relative to one another. 
10 - 
In order to produce an image in which the details of the 
article are shown, rather than the pure distance 
information, the processed image can then be applied to 
the intensity image taken at plane I 0 which is the in focus 
15 plane, to give both the intensity information of the image 
and the depth information of the image in the one image. 
Thus, the intensity information provides the structural 
detail of the various components in the image and the 
depth information, as previously described, enables you to 
ascertain the relative position of parts of the image. 
Thus # this image will not only contain the depth 
information, but will also make it much easier to 
ascertain what various components of the image are, 
because the actual structure in the image will be seen. 

Thus, according to the preferred embodiment of the 
invention, the imaging range of the lens on the camera 20 
and the colour rendering can be calibrated and therefore 
positions of the image can be quantified. This technique 
then offers the opportunity to provide a passive ranging 
capability. in conditions of near thermal uniformity, 
where there is very little contrast in the field of view, 
the depth positional information provided by this process 
suggests another contrast mechanism for viewing and 
processing images of scenes. Finally, under obscured 
viewing conditions, such as under camouflage, the outline 
of an object can be disrupted, making identification 
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difficult. The colour depth rendering of the scene allows 
the piecing together of image components at common 
positions, as even though the outlines can be disguised, 
the physical positions in space cannot be* 

The manner in which the algorithm operates to process the 
variance data to produce the phase image is described 
below. 

At each point in space, an optical beam possesses two 
properties: intensity and phase. Intensity is a measure 
of the amount of energy flowing through each point r while 
phase gives a measure of the direction of the energy flow. 

Intensity may be measured directly, for example by 
recording an image on film. Phase is typically measured 
using interference with a "reference beam". in contrast 
the present method gives a non- interferons trie method for 
measuring phase. 

Intensity can be measured over two parallel planes V + 
extending across the direction of propagation of the wave 
field on the side remote from the incident radiation. 

The present invention determines phase by providing a 
solution to the transport-of -intensity eolation: 

(1) V x .(/V x <* ) = 

where I is the intensity in the plane, the gradient 
operator in the plane is denoted V ± , fc is the wave number 
of the radiation, and dl/dz is the intensity derivative or 
rate of change of intensity. However, in order to obtain 
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the depth information which is previously referred to, 
rather than using the intensity values which are obtained 
by the camera at the two planes, the variance values 
calculated in the above mentioned manner are obtained so 
that a variance derivative, or rate of change of variance 
dVIdz. Note that dV/dz is estimated from the difference of 
the measurements in the planes V. and v + shown in Figure 1, 

while the intensity X is given by the average of the 

measurements. 

In order to obtain a solution to equation 1 the function A 
is first defined as: 

(2) V ± A = W ± 0 . 

Thus equation (1) becomes: 

» 

(3) V ± .(V x A)=-kd z V . 

Making use of a standard identity V ± • V ± = V x 2 , this may be 
written: 

(4) V ± 2 A = -kd z V 

where V ± 2 denotes the two-dimensional Laplacian acting 
over the surface of the image. This equation has the 
following symbolic solution: 

(5) A = -kV x ~ 2 d z V. 

If the gradient operator V ± is applied to both sides of 
this equation, it becomes: 
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(6) V ± A = -kV x V x '' 2 d z V . 

The defining equation (2) for the function A allows (6) to 
be transformed into: 

(7) IV ± 0 = -kV ± V L - 2 d z V . 
Dividing both sides by I then yields: 

(8) v A #=^f-«v J .v;- a a f v. 

However, as is mentioned above, rather than dividing both 
sides by the intensity value X, a maximum focus value 
across all of the images could be used, or the variance 
data at the plane could be used, or a maximum variance 
across all of the variance state of values could be used. 

Taking the two dimensional divergence V x • of both sides of 
(8), and again making use of the standard identity 
V ± «V X =V X 2 , then (8) becomes: 

( 9 ) V x > = -fcV x • J/" 1 V x V X ~ 2 3^J . 
This equation has the following symbolic solution: 

do) ^-^v x (v x 4rv x v x - 2 3 z yj). 

In order to implement a practical solution to equation 
(10), the following formulae are required • A suitably- 
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well-behaved function f (x,y) may be written in the form of 
a two-dimensional Fourier integral: 

(11) /(*,)-] J hk x ,k y )e iik <** k > y) dk x dk r 

— OO —OO 

A 

The function f{k x ,k y ) is called the "Fourier transform" of 

f(x,y). 

The x derivative of (11) yields s 



< 12 > J J [ *«/( ^r.*, )V k ^ y) dk x dk r 



Hence the Fourier transform of — f (x,y) is equal to the 

Fourier transform of f(x,y) multiplied by fit,. Stated 
d 

differently, — = /F" 1 fc x F # where F denotes Fourier 

ox 

transformation and F" 1 denotes inverse Fourier 
transformation* Similar considerations apply to 

d 

-r-f (x # y) • 

If the Laplacian V 2 ± = —+— of (11) is obtained and 

ox oy 

similar reasoning applied, it follows that VT 2 =-F" 1 & ~ 2 F 

J. f f 

where . k? =k x l +k y 1 . Thus : 



(13) V?=-F-VF, = ~ = iF-%P 

dx 3y y 
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Here, F denotes Fourier transformation, F~ x denotes inverse 
Fourier transformation, (k* fcy) are the Fourier variables 
conjugate to (x,y), and 



Equations (13) can be used to rewrite equation (10) in the 
form 



In practice division by intensity is only performed if 
that intensity is greater than a certain threshold value 
(eg. 0.1% of the maximum value). 

Division by k r does not take place at the point kr = 0 of 
Fourier space; 

instead multiplication by zero takes place at this point. 
This amounts to taking the Cauchy principal value of the 
integral operator V ± ~ 2 

In order to quantitatively measure the phase of object it 
is necessary to incorporate some physical constants into 
the phase recovery algorithm given in Equation (14) 
relating to the experimental setup in use to quantify the 
variables k** ky. This can be done by rewriting equation 



5 



kf~k* +k* . 



10 



(14) 
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(14) in the following form suitable for implementation 
using a fast Fourier transform: 



lit 1 



F 1 



A&(NAxf i 2 + f /(x,y) 



i 2 + f 



2k 



^(iVAx) 2 i 2 +/ /(x,y) F + j 



where i 



. \-N N~\ . 



index the frequent components of F(V + — VL) 



where the intensity derivative d z V(x,y)L& obtained by 
subtracting two images V+ and VL separated by a distance 
Szt i and j are the pixel numbers on the image, and using 
the fact that the Fourier space step size is given by 



Ak = - 



NAx 



where the image is an N x N array of pixels of size Ax . 
Thus in addition to measuring the three intensity 
distributions it is necessary to know the pixel size Ax , 
defocus distance Sz and wavelength X in order to make a 
quantitative phase measurement. All of these quantities 
can be readily determined: the pixel size can be 
determined directly for example from the CCD detector 
geometry (in the case of direct imaging) , or by existing 
techniques for calibrating the transverse distance scales 
(in the case of an imaging system) , the defocus distance 
can be measured directly, and the spectral distribution of 
the illumination can be determined either by 
monochromating the incident field or by analysing the 
spectral distribution of the radiation using existing 
spectroscopic methods* 
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An example of the phase-retrieval method implementing the 
solution of equation (14) can be represented by the 
flowchart shown in Figure 6. As shown in Figure 2 the 
quantitative determination of phase of a radiation wave 
field proceeds from a set of variance measurements {V B } 
over the two spaced apart planes V_ and V + . A measurement 
of central intensity I{x # y) in a selected plane parallel 
to and midway between the planes V_ and V + is also 
obtained. The intensity measurements are performed over a 
defined array on each of the two planes A and B and the 
respective values subtracted to produce a measure of the 
intensity derivative. This value is multiplied by the 
negative of the average wave number. The data are split 
into two component sets and a fast Fourier transform is 
performed to produce the respective x: and y components in 
the Fourier domain. A filter is then applied to the 
Fourier domain representations to correspond to the 
inversion of a differential operator reflected in the 
untrans formed representation. An inverse Fourier 
transform is then performed o*i each of the x and y 
components to produce a spatial domain value from which 
the differential operator has been removed. A division by 
the central intensity I(x,y) obtained by averaging the 
intensity measurements over planes V + and V- is then 
performed if the intensity level is above a selected 
threshold level. The resultant data is again Fourier 
transformed and multiplied by the same filter to again 
correspond to the inversion of a differential operator 
reflected in the u,ntrans formed data. The resultant 
components are again inverse Fourier transformed and 
summed to provide a retrieved phase measurement . 
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It will be apparent that in general the method according 
to this invention can proceed from any suitable 
representative determination of intensity derivative or 
rate of change of intensity over a selected surface 
extending across the propagation direction and the 
intensity over that same surface. As will be explained in 
various examples these data can be obtained in a variety 
of ways and the method implemented to yield phase of the 
radiation wave field. 

Rewriting equation (14) with: 

&Ak x ,k y9 a)~k x k r ~ 2 
£l y (k x ,k y ,a) = k y k~ 2 

Mx;y)=0 M (x,y)+0 W (x,y), 

gives 

(15) 

0 <x) (*,y)=F-^ 

where : 

<p(x,y) denotes the recovered phase, 

F denotes Fourier transformation, and F' 1 denotes 
inverse Fourier transformation, 



H:\Karen\Keep\Iatia Prov - Image containing depth informationl.doc 27/02/04 



X(x,y) is the intensity distribution over the plane 
of interest, 

(x,y) are Cartesian coordinates over the plane of 
interest, 

(k x7 k y ) are the Fourier variables conjugate to (x,y) 

k = 27tj X is the average wave number of the 
radiation, 

/L is the average wavelength of the radiation, 
dV/dz is the estimate for the longitudinal variance 

derivative, 

a is the regularization parameter used to 

stabilize the algorithm when noise is present. 

As given above, the solution to the transport of intensity 
equation (1) assumes a perfect imaging system. That is, 
there are no ^aberrations" present in the optical system 
used to obtain the intensity data which is fed into the 
algorithm. Of course, no imaging system is perfect. The 
imperfections present in an imaging system may be 
quantified by a set of numbers: 

(10) A l ,A 2 ,A 39 .>. 

which are termed aberration coefficients* 

If intensity data were taken on an imperfect instrument 
whose imperfections were characterized by a certain set of 
known aberration coefficients A lJ A 2J A^^. , it would be 
desirable if the filters £l x (k x k y9 a) and Q, y (k x k y ,a) present in 

(15) could be replaced by modified filters which 
explicitly depend upon the aberration coefficients: 
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(17) & x (* x k r ,a,A t ,A 2 ,A 3 ,..~.) and ^(M^A^A,,...) 
This would allow the imperfections of the imaging system 



quantitatively correct phase retrieval using aberrated 
imaging systems. For the special case of a non-absorbing 
phase object in a radiation wave field of uniform 
intensity with weak (i.e. much less than In radians) phase 
10 variations the appropriate modified filters lead to the 
following functional form for the phase-retrieval 
algorithm: 



5 



to be explicitly taken into account, leading to 



(18) 



15 



where : 



Vkberratea (x# y) is the aberrated variance measured at 
defocus distance Sz , 

Ann are the aberration coefficients which 
characterize the imperfect imaging system. 



20 



If a filter is defined: 



.(19) &(k x k y ,a,A l ,A 2 ,A 3 ,...) 



1 



25 



Then (18) becomes: 



(20) 



Hx, y) = F-^F-iF'^Fi^ (x, y) -1 } 



«o 
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The term { I^^^y)-! }is a measure of rate of change of 
intensity. x 0 intensity is a measurable constant for 
uniform intensity so that (20) is the same general form as 
(15). Consequently the special case of aberration can be 
i dealt with by changing the filter in the general method 

described above. The x and y component filters £l x and Q y 

are given by 

(21) Q,=Q,=lfl 

V2 

Figure 7 schematically shows an arrangement for 
quantitative phase amplitude microscopy. A sample is 
illuminated using a source of white light Kohler 
illumination 15, commonly found on optical microscopes. 
The light is transmitted through an object 16 and 
collected by a microscope imaging system 17 and relayed to 
a CCD camera 18 or other digital imaging device having a 
planar imaging surface. Three images are collected: an 
in-focus image, l c , and two slightly out of focus images I + 
and I_. The defocus is obtained by suitable means such as a 
drive system 19 to adjust the microscope focus knob. The 
defocus introduced is usually quite small. so that 
degradation in spatial resolution is minimised, although 
the optimal amount of defocus to use is determined by 
sample properties and imaging geometry such as 
magnification, numerical apertures, etc. 

When taking the images the numerical aperture of the 
condenser is chosen to be less than the numerical aperture 
of the objective being used. if this is not the case then 
serious image degradation will occur, although the precise 
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amount by which the condenser and objective numerical 
apertures should differ involves a tradeoff between image 
fidelity and spatial resolution with the optimal 
difference depending- on the sample properties and the 
optics used. 

The variance data values are determined in the manner 
previously described from the intensity values collected 
at the planes V + and V. and are subtracted to produce a 
representative measure of rate of change of variance 
(variance derivative) . From this and the intensity data 
of collected image I 0 the method described above can be 
used to produce quantitative information about the 
magnitude of the phase shift in the image plane. 

There may be cases in which it is desirable to take more 
than two images in order to obtain a better estimate of 
the intensity derivative dV/dz. A function can then be 
fitted to this data from which dV/dz can be computed and 
used in the phase determination method in place of the 
simple subtraction of two images normally used. 

It is also possible to operate this system in reflection 
geometry to obtain surface topography. The principle of 
operation is the same, however the optics have to be 
folded back on themselves to form a reflection geometry - 
otherwise the process is identical. 

For certain applications it can also be desirable to 
filter the light to a particular wavelength, although this 
is not necessary for the described imaging process as it 
works equally well with white light. 



H:\Karen\Keep\Iatia Prov - Image containing depth informationl.doc 27/02/04 



- 31 - 

An implementation is shown in Figure 8. An Olympus BX-60 
optical microscope 20 was equipped with a set of UMPlan 
metallurgical objectives and a universal condenser to 
provide Kohler illumination. In order to be able to 
compare the results with existing imaging modes Nomarski 
DIC optics and a set of cover-slip corrected UplanApo 
objectives were also acquired for this microscope, 
enabling images to be taken of the same field of view 
using both phase retrieval and Nomarski DIC for the 
purposes of qualitative comparison* A 12 -bit scientific 
grade Photometries SenSys CCD camera 21 equipped with a 
1300x1035 pixel Kodak KAF-1400 CCD chip was added to the 
0.5x video port on the microscope to acquire digital 
images of the sample • 

The phase recovery technique of this embodiment of the 
invention requires the acquisition of defocused images. A 
stepper motor drive system 22 was attached to the focus 
knob of the microscope. This stepper motor 22 was 
attached to the parallel port of a 133MHz Pentium PC 23 
also used to control the CCD camera 21, enabling full 
automation of the acquisition of through- focus image 
sequences. This data acquisition system was linked to 
custom software written to recover phase images from the 
CCD images, thereby enabling full automation of the image 
acquisition and data processing sequences. 

The form of the differential operators used in the 
preferred embodiment of this invention are based on the 
optics of the system used to obtain the above-mentioned 
images. Thus, the operator takes into account the details 
of the optical system used to take the images. This is 
achieved by: 



H:\Karan\Keep\Iabia Prov - Image containing depth informationl.doc 27/02/04 



Determine the numerical aperture of the condenser, NAcondcnsor 

Determine NA of objective, NAo b j ective # and pobjective, the radius 
of the objective aperture 



NA , 

z= condenser 



(These are settings and dimensions on the microscope.) 

Determine radial lateral spatial frequency, p, and 
longitudinal spatial frequency, T|. 

(These are dependent on the pixelation and position 
distribution of images taken in the series.) 

Determine the wavelength, X, of the radiation to be used. 
The new form of the operator is 



T p +a 2 
where, 

T P (p) = 2m Sz faT™ (p, 7j)d7j 
and 



1* 



2np 



l -[ 1 ^^ +1 )-^-(^J-(l + ^(f-lfl 



Figure 9 is a flow chart generally illustrating how T p It 
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determined by means of the above equation merely showing* 
breakdown of the various components of the equation. 

The calculator depth image can then be presented in 
different modalities to give better visualisation. Such 
modalities including DIC, bright field images, dark field 
images and other conventional modality images. Techniques 
for forming these types of images from the phase data 
determined in accordance with the present invention are 
described in our copending International Application No. 
PCT/AU02/00590 (the contents of which is incorporated into 
this specification by this reference) . 

Since modifications within the spirit and scope of the 
invention may readily be effected by persons skilled 
within the art, it is to be understood that this invention 
is not limited to the particular embodiments described by 
way of example hereinabove. 

In the claims which follow and in the preceding 
description of the invention, except where the context 
requires otherwise due to express language or necessary 
implication, the word "comprise", or variations such as 
"comprises" or "comprising", is used in an inclusive 
sense, ie. to specify the presence of the stated features 
but not to preclude the presence or addition of further 
features in various embodiments of the invention. 

Dated this 27th day of February 2004 

IATIA IMAGING PTY LTD 

By their Patent Attorneys: 

GRIFFITH HACK 

Fellows Institute of Patent and 
Trade Mark Attorneys of Australia. 
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